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Abstract 

Two discretizations of the vector nonlinear Schrodinger (NLS) equation are studied. One of these 
discretizations, referred to as the symmetric system, is a natural vector extension of the scalar integrable 
discrete NLS equation. The other discretization, referred to as the asymmetric system, has an associated 
linear scattering pair. General formulae for soliton solutions of the asymmetric system are presented. 
Formulae for a constrained class of solutions of the symmetric system may be obtained. Numerical 
studies support the hypothesis that the symmetric system has general soliton solutions. 

PACS Numbers: 03.40.Kf, 46.10+z, 42.81.Dp, 42.81. Gs 

1 Introduction 

In recent years there has been wide interest in the study of solitons and integrable systems. Researchers 
have found that not only are continuous systems (i.e. PDE's) integrable via the inverse scattering transform 
(1ST) but also that interesting classes of discrete systems (semi-discrete as well as partial difference equation) 
are integrable - cf. Ablowitz and Segur (l) for an early review). In this letter, we discuss soliton solutions 
and the integrable nature of certain discrete systems associated with the vector extensions of the nonlinear 
Schrodinger (NLS) equation: 



iq t = Qxx + 2q\q\ 2 . (1) 

NLS is a centrally important and physically significant nonlinear equation which possesses solitons and is 
integrable via 1ST ||. Furthermore, NLS arises in many areas of physics, such as the evolution of small 
amplitude slowly varying wave packets in: deep water, nonlinear optics and plasma physics (see e.g. Q). 
In 1974, Manakov f| showed that the vector NLS (VNLS), 

iq_t = <\xx + 2 ||q|| 2 q (2) 

where q is an A-component vector and ||-|| denotes the vector norm, also possessed solitons and could be 
integrated via 1ST (actually, in || only the case N = 2 was studied in detail; however the extension to the 
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N-th order vector system is straightforward). The second order (N = 2) VNLS equation is relevant in the 
study of electromagnetic waves in optical media in which the electric field has two nontrivial components. In 
optical fibers, the components of q in eq. (|^) correspond to components of the electric field transverse to the 
direction of wave propagation. These components of the transverse field compose a basis of the polarization 
states. Although, in optical fibers, the equations governing the field are in general a non-integrable variation 
of vector NLS [Q , there are circumstances in which eq. (|^) is the appropriate model || ^) . 

In 1976, Ablowitz and Ladik |Q found that the following discrete system is integrable via 1ST: 



i — In = T7(9«-i - 2q n + q n+1 ) + \q n \ 2 {q n +l + q n -l) (3) 
At h z 

The continuum (h — » 0) limit of this discrete system (||) is NLS (|]) hence this system can be referred to as 
integrable discrete NLS (IDNLS). The discrete system itself is also useful in physical applications (see e.g. 
111,0). 

IDNLS provides an excellent numerical scheme to solve NLS (|l|) (cf. fl2[| ). Computations demonstrate 
that, by simply replacing the nonlinear term in ([|) with 2\q n \ 2 q n , one obtains a poor numerical scheme. 
In contrast, the preservation of the integrable structure in IDNLS (J3J) plays a key role in its utility as a 
numerical scheme. This points out that care must be taken in the choice of discretization of the nonlinear 
terms for the purpose of numerical simulations. 

It is natural to look for useful discretizations of VNLS (^). However, the special character of IDNLS 
as compared to other discretizations of NLS, suggests that different discretizations of VNLS may have very 
different dynamics. Therefore, the choice of discretization merits close analysis. 

In this letter, we discuss two discretizations of VNLS which we refer to as the symmetric and asymmetric 
discretizations. These discretizations are as follows: 

Symmetric discretization: 



d 1 2 

i— q« = 7j(q«-i - 2q„ + q„+i) + q n |] (qn-i+qn+i) (4) 

dt ll z 



where q„ is an ^-component vector. 
Asymmetric discretization: 



i-jq™ = To (q«-i - 2q„ + q„+i) - (r£q n _i)q„ - (r^q„)q„+i (5a) 



d _ 1 

where, as before, q„ is an TV-component vector as is r„. The superscript T denotes the transpose. 



-i^jr n = 72 (r n _i - 2r„ + r„ +1 ) - (r^q„)r„_i - (r^ +1 q n )r„ (5b) 



Both systems (Q, |5a|f3b| ) reduce to VNLS (0) in the continuum limit (h — ► 0) where in (|5a[-pb|) we take 
r = — q*, where * denotes the complex conjugate, after taking the continuum limit. The asymmetric system 
(|5a ■ 5b ) , however, does not admit the symmetry r„ = — q* and hence remains "asymmetric" for h finite. 

These systems have a number of interesting properties. Both systems (^, [5a-5b) reduce to the integrable 
discrete scalar case (||), just as VNLS (|^) redu ces to NLS (0). The asymmetric system ([5a|-5b) is associated 



with a linear operator pair (see eq. lla- llb| ). Although there is no known linear operator pair for the 



symmetric system (^), under the reduction q„ = e M 'v n , where v„ is independent oft, this equation reduces 
to an V-dimensional difference equation which is known to be integrable (cf. |]l3| ). 

Both the symmetric (^) and asymmetric (|5a|-pb|) systems possess a class of soliton solutions. The asym- 
metric system has a class of soliton solutions which correspond to all the soliton solutions of the continuum 
limit NLS (J2J) . However, the exact multi-soliton solutions found so far for the symmetric system (|~ 
reduce to a subset of those associated with eq. @ . Since the known exact soliton solutions of eq. 
only a subset of those known in the continuum limit (j^) , we numerically examined the interactions of solitary 
waves associated with the symmetric system which lie outside this class. The numerical evidence indicates 
that, in the symmetric system, the solitary waves interact elastically and are therefore true solitons. This 
finding supports the possibility that the symmetric system (0) is indeed integrable. 
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2 NLS, Integrable discrete NLS and Vector NLS 



We briefly review the relationship between NLS, IDNLS and Vector NLS in order to clarify their relation to 
the symmetric and asymmetric discretizations of vector NLS. 

2.1 NLS and IDNLS 

NLS (|l|) can be obtained as compatibility condition of a linear operator pair (e.g. a Lax pair). The linear pair 
can be used to solve (0) by 1ST. The existence of an associated linear operator pair and 1ST are hallmarks of 
integrability. Therefore, in the effort to construct an integrable discretization of NLS, it is natural to look for 
a discrete version of the operator pair and the compatibility condition. In fact, the integrable discretization 
of NLS (||) was obtained in exactly this manner |jj . 

The extent to which the integrable structure of NLS is preserved under the integrable discretization is 
significant: The 1ST for IDNLS is analogous to the theory for NLS. In particular, the one-soliton solution of 
IDNLS (|) given by 

q n (t) = Ae- l{bhn -^-^sec\i{ahn -vt-6) = e*q n (t) (6) 

with 

smh(ah) 2(1 — cosh(ah) coslbh)) s'mh(ah) sixiCbh) 
A= ^^' U = J? ' V = 2 h* 

where a,b,8,4> are arbitrary parameters, converges (as h — > 0) to the one-soliton solution of NLS with the 
same four free parameters. More generally, solutions of IDNLS converge to the solutions of NLS with an 
error of 0(h 2 ). 

For both NLS and IDNLS, generic rapidly-decaying initial data resolves itself into some number solitons 
plus radiation that vanishes as t — > oo (in the sup norm). Even though the solitons interact nonlinear ly 
upon collision, it is well-known that they retain their shape and speed after the collision. This can be 
determined by comparing the forward (t — > +oo) and backward (t — > — oo) long-time limits of the solutions. 
In these long-time limits, solitons traveling at different speeds are well-separated. Although this collision is 
elastic, the position and overall complex phase of the individual solitons are shifted between the forward and 
backwards long-time limits. This elastic interaction of solitons up to a shift in phase (position and overall 
complex phase) is typical of (1 + l)-dimcnsional integrable systems. 

2.2 VNLS 

VNLS (||) can be obtained by substituting the vector q for q in NLS ([!]) and replacing complex conjugation 
with the Hermitian conjugate. Similarly, the linear operator pair for VNLS can be obtained by appropriately 
making the same substitution of scalars by vectors in the linear operator pair for NLS (cf . Manakov || ) . 

When considered as a model for the propagation of electromagnetic waves in optical fibers, the compo- 
nents of VNLS play the role of a basis for the polarization vector. In the derivation of VNLS by Manakov || 
the choice of basis is arbitrary. Therefore, the vector system (|^) ought to be, and is, invariant under a change 
of basis for the polarization. Mathematically, a change of basis is obtained by multiplying the independent 
variable q by a unitary matrix, U. Therefore, the freedom in the choice of basis is reflected in the fact that 
the vector system (|^) is invariant under the transformation q — > Uq. This symmetry is an important feature 
in distinguishing among discrete versions of vector NLS (notably, both the versions of discrete vector NLS 
discussed in this letter fl [5a|]5b|) retain this symmetry). 

Under the reduction 



Q = cq (7) 

where c is a constant, TV-component vector such that ||c|| 2 = 1 and q is a scalar function of x and t, VNLS 
reduces to NLS (pi). This is a manifestation of the fact that the vector system is a generalization of 
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NLS obtained by allowing the polarization to be non-constant. When the solution of VNLS has a constant 
polarization, NLS is recovered. As a consequence, any solution of scalar NLS has a corresponding family of 
solutions of VNLS. We call a solution of the form (Q) a reduction solution. Because the vector c is arbitrary 
(up to the constraint that it is of unit length), a reduction solution has correspondingly more parameters 
than the solution of scalar NLS to which it reduces. We refer to the vector c in eq. (Q) as the polarization 
of the reduction solution. 

The soliton solutions of vector NLS can be found by 1ST || . The one-soliton solution of vector NLS is 
the reduction solution (i.e of the form (Q) where q n is the one-soliton solution of NLS. As in the scalar case, 
the 1ST H shows that solitons interact elastically. By elastic, it is meant that the solitons retain their shape 
and speed after interaction as shown by comparing the forward and backward long-time limits. 

For the vector system, care is needed: the nature of the phase shift for vector solitons is somewhat 
different than the scalar case. In addition to a shift in the center of the peak (as in the scalar case), vector 
solitons also undergo a change in polarization upon collision. That is, in the forward and backward long-time 
limit, the solution asymptotically approaches a linear superposition of the individual solitons: 



q ~ 2_, 'if as i-> ±00 

j 

where q^ = cfqf and, for each j, q^ and qj are one-soliton solutions of NLS with the same amplitude and 
speed. Comparison of the forward (+) and backward (— ) long-time limits shows that 



but 



||c+|| 2 H|c7|| 2 = l. (8) 

A closed formula for this shift in the polarization can be calculated by considering the eigenfunction of the 
associated scattering problem Q . 

We refer to the squared modulus of a component of the polarization vector as the intensity of that 
component of the polarization. That is, if Cj is the £-th component of Cj, the polarization vector of the 

soliton j, then the intensity of the of the £-th component of soliton j is |c^| 2 . The relation (||) implies 
that the total intensity of each soliton is preserved, specifically the sum of the intensities is equal to one. 
However, the distribution of intensity between the components of the polarization of an individual soliton 
will not, in general, be equal in the forward and backwards long-time limits due to interaction with other 
solitons: subject to the constraint of eq. (||), in general, 



This change in the distribution of intensity is a distinctive feature of the vector system. There is no cor- 
responding phenomenon in the soliton interactions of the scalar equation (in both the scalar and vector 
equations, the location of the peak of an individual soliton is shifted by the soliton interaction). 

Manakov H also observed the following special case for the interaction of the solitons of VNLS: If, for 
every pair of solitons j and k, either 



(9a) 



(9b) 
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where • denotes the dot product, then 



(10) 



for all j. As a consequence of (10), \cp \ 2 = \c^'^\ 2 for all j,£. Physically, this means that, in any 



-W+I2 



polarization basis, the distribution of intensity among the components of polarization for an individual 
soliton is not changed by the soliton interaction. 



3 Asymmetric Discrete Vector NLS 

The preceding derivations of IDNLS and VNLS from NLS suggest two methods to obtain an operator pair 
which has a discrete form of VNLS as its compatibility condition: (i) discretize the operator pair for vector 
NLS in a manner analogous to that used to obtain the operator pair of IDNLS or (ii) appropriately substitute 
vectors in to the operator pair for IDNLS. These approaches yield the linear operator pair 



zl N hq n 



T ■ 1 ( 

! 1'> r n-l — l 2h 2 \ Z ~ 



ZTt 



-4 (Jc 

-ir^q„_i 



.. - zq n -ij 



(lla) 
(lib) 



where I at is the N x N identity matrix. The asymmetric system (|5a-5b) is the compatibility condition of 
this pair. 

Under the reduction 



q„ = cq n , r„ = dq* (12) 

where c and d are constant vectors such that c • d = 1 and q n is a scalar function of n and t, the asymmetric 
system (|5a-5b) reduces to IDNLS (||). Solutions of this form ( p"2| ) are referred to as reduction solutions. 



Hence, every solution of IDNLS generates a family of reduction solutions of the asymmetric system. We are 
particularly interested in reduction solutions where d = — c* because, in this case, r„ = — q* and (§a|-ft]) 
reduces to a single equation. 

3.1 General Formula for Multi-soliton Solutions 

Multi-soliton solutions of the asymmetric system can be derived, without the use of the 1ST machinery, by 
Hirota's method. First, make the independent-variable transformations 

Sn Sn 

qn — , r n — — 

In Jn 

where g„ and g„ are vectors with N components and f n is a scalar. Then, solutions of the system of bilinear 
equations 



/n+l/n-l ~fn = SnSn (13a) 

ih D t fn ■ gri = /n-lgn+1 ~ 2/„g„ + /ri+lgri-1 - H„g„ (13b) 

-ih 2 D t f„ ■ g n = /n-lgn+1 ~ 2 /„g„ + /n-lgn+1 - H„g„ (13c) 

/„H„+i = h 2 (g n gl +1 - Sn+lgn) (13d) 

/„H n -l = h 2 (gn-lSn - SnSn-l) ( 13e ) 

where the N xN matrices H„ and H„ are auxiliary variables, are solutions of the asymmetric system (pi -5b). 
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The M-soliton solutions for these bilinear equations (13a- 



13e) are given by the following determinants: 
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where: I is the M x M identity matrix; 

(i) A„ and B are M x M matrices defined by 



,h(pt+P* j ) 



with 



>li,, 



: Pj nh +j^(2-e 



- h P i 



(14a) 



(14b) 



(15) 



the complex numbers pj = aj — ibj, aj > 0, determine the amplitude and speed of the j-th soliton; the 
iV-component complex vectors jj and 'y- (J = 1, ... , M) determine the polarizations and envelope phases 
of the solitons; 

(ii) $ n are, respectively, the M-component column and row vectors 



$W = e -hPe e m,n <&G) = e v], n . 
(iii) are, respectively, the M-component row and column vectors, 



(j) 



3.2 The One-Soliton Solution 



"7- 



(fc) 



-hp* e -(k) 



From the above formulae (14a-14t), the one-soliton (M=l) solution is 



_^ e _ ihbl smluVO e ^,„ sech(ai + ^ 



r„ = — 



(7?7i) 5 ft 



where 



(16a) 
(16b) 



0i <n = Im »7 1)n 



sinh(ai/i) sinf&i/i) 
ai/m — 2 — i 



-b\hn 



h 2 

2(1 — cosh(ai/i) cos(6i/i)) 
h 2 



(17a) 
(17b) 
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and 



£1 = log 2 {tT 7i} - log i 1 — ^ 



Note that to get a solution (16a- 16b| ) with 5\ real we restrict our attention to the case where 7^7! is real 
and positive. The one-soliton solution (16a -16b) is a reduction solution- i. e. of the form ( |l2] ) -where 



7i 



-ihb\ 



7i 



ihb\ 



(7f7r) 



and the scalar function q n is given by eq. (g)- the one-soliton solution of IDNLS -where a = a x , 6 = b x and 
= —Si. In this solution, d = c* if, and only if, j x =7*. 



3.3 Two soliton Interaction 



To show that, for M > 1, the determinants (14a- 14bj ) indeed give a multi-soliton solution, we consider the 
long-time limits. In these limits, solitons moving at different speeds are separated. For concreteness, let 
M = 2 and 



sinh(aift,) sin(6i/i) > sinh(a2/i) sin^h) 



(18) 



for the given 01,61,02,62. The condition ( |l8| ) ensures that each soliton travels with a different speed (the 
analysis below holds, in general for M solitons as long as each has a different speed). When the solitons 
travel at different speeds, we determine the asymptotic form of an individual soliton by taking long-time 
limits in a coordinate frame moving with that soliton. 
The limits 



sinh(a,/i) sin(6,/i) , . 

n ~ 2 y J \ v 1 J t, t^±oo (19) 

CLjh 4 

are long-time limits in a coordinate frame moving with soliton j. In the limit ([l^) with j = 1 



Re r]i in = const, and Re r\2. n ~ * ±00 as t — > ±00 



(20) 



The substitution of (|20|) into eq. (14a- 14b), where M — 2, yields 



and 



as t — > ±00, 



(21) 



which is the asymptotic form of soliton 1 in the forward (+) and backward (— ) long-time limits. The long- 
time limits for soliton 2 are similar. The coordinate frame of soliton 2 is obtained by the the limit (|l9|) with 
j = 2 which yields 



Re 371, 



^00 and Re 772.-, 



const 



as t — > ±00 



(22) 



Note the change in relative sign between 77 l ra and t in the coordinate frame of soliton 2- the slower soliton 
-in the long-time limit. To obtain 



and 



L 2,n 



as t — > ±00, 



(23) 
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which is the asymptotic form of soliton 2 in the forward (+) and backward (— ) long-time limits, substitute 
(p2|) into (|14a|-|14b|), where M = 2. Combining (J2lj) and (p3|) gives 



The effect of interaction on the solitons is determined by comparing q~ n and r~ n with and r+ n for 
j = 1, 2. We give the formulae for the asymptotic forms of the solitons below. 
First, we consider soliton 1. In the backward long-time limit, 



1l,n — C l 9l,n' 



l,n — ^1 #l,n 



where 



c, = 



7i 



(7f7i) 



-ihb\ 



(7i7i) 



and n is a one-soliton solution of IDNLS (the <j in eq. (§)) with a = <zi, b = b\ and 



(24a) 
(24b) 



= -*r = - lo s 2 {7^71} + log I ^ — 



(25) 



As in the one-soliton case, we require that is real and positive. The forward long-time for soliton 1 is 



(26) 



where 



-(e~ h ^a 2 - p*) (a 2 (7^7 2 )7i - P*(^li)~f2 

X 

1 

e 

X 



-hpl 



(cr 2 - p) (^2(72 7 2 )7i - ^(7271)72) 



and qi n is a one soliton solution of IDNLS with a = ai, b = b\ and 



= -<^ + = - log'X + logcr 2 + log {1I12} 



also 



(27a) 
(27b) 



X = {(<?i<?2 |p| 2 ) (a 1 a 2 (7f7i)(7r7 2 ) - \p\ 2 (lIl2)(lIli))V 

and we require that 7j,7j f° r J = 1)2 are such that x anc l 72*72 are rea, l anc l positive. The constants cr\, 
cr 2 , p are given by 



ci 



z Hpi+pi) _ ^ : 



(T 2 



e h(p2+Pi) _ ]_ : 
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Hence, the expression {o\<Jri — |p| 2 ) is real and positive. Note that, although c 1 ^ c^, d 1 ^ d^, and 
6± ^ 81 the parameters a\ and b\ are the same in both long-time limits. Thus, soliton 1 has the same form 
in both long-time limits but undergoes a phase shift due to interaction with soliton 2. 

The asymptotic limits of soliton 2 are similar to those for soliton 1 except that the relative change in 
sign in (E2n reverses the calculations, that is, 



c 2 ?2, 



and 



h 9 2 ,-. 



where c 2 has the form (27a) and d 2 has the form (27b) with the indices 1 and 2 exchanged. The scalar 



function q 2 n is a one-soliton solution of IDNLS with a = a,2, b = ri2 
with the indices 1 and 2 exchanged. Similarly, 



and 



= — S 2 where S 2 is equal to 5i 



l2,n — C 2~<Sn 



and 



L 2,n 



where c 2 has the form (24a) and d 2 has the form ( p4b| ) with the indices 1 and 2 exchanged and q 2n is a 
one-soliton solution of IDNLS with a = 0,2, b — ri2 and 9 = —8 2 where 8 2 is equal to J-f with the indices 1 
and 2 exchanged. 

In order to construct a solution such that 



r j,n %',n 

for j = 1,2 it is necessary and sufficient that 



7i = 7i 



7 2 = 



'iTi^l 2 + |e-' lai 7l 2) | 2 (1 - e- 2 '^)7 1 1) * 7 f ^ 



hif (1 - e-^) 7 W 7 f* 



e -hai 7 (l) 



7i 



(2; 



75- 



(28a) 
(28b) 



For this choice, x an( i 7j7j f° r J = 1>2 are real and positive as we required . However, if ( |28a -28b) hold 
then 



4n * -C ( 29 ) 



for j = 1, 2. In fact, ||r+ n — (— q+*)|| = O(h). An alternative is to set 



7i=7* and 7 2 =72- 



(30) 



Then, if (30) holds, x and jjlj- for j = 1, 2 -are real and positive and 



— 42, n 



but 



K„ - (-q^)ll = o(h) = i|r 2 -„ - (-q 2 -;)||. 

In contrast, recall that the solutions of the scalar IDNLS converge to solutions of the scalar NLS with 
0(h 2 ) convergence. This unavoidable asymmetry of the system (j^[5b|) makes it less desirable as a discrete 
approximation of VNLS. Hence, we turn our attention to the symmetric system. 
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4 Symmetric Discrete Vector NLS 



The symmetric system (Q) is the natural vector gener alization of IDNLS (§). We now describe some important 
symmetries of the symmetric system. In analogy with the scalar case (^), the symmetric system can be 
thought of as the reduction of the system 



i—;<in = To (q™-i - 2 q™ + q«+i) - r «qn (q«-i + q«+i) (3ia) 

-i^r„ = (r„_i - 2r„ + r n+1 ) - q^r„ (r n _i + r n+1 ) (31b) 

obtained by letting r„ = — q* , which is the symmetry that is broken by the asymmetric system. 
Under the reduction 



q« = cq n (32) 

where ||c|| = 1, the symmetric system reduces to the scalar IDNLS (||). Furthermore, if q„ satisfies the 
symmetric system, then so does Uq„ where U is a unitary matrix (this symmetry is a discrete form of 
eq. (0)). Therefore, the symmetric discretization retains the reductions and symmetries of the asymmetric 
system and has the additional symmetry r n = — q*. 

Despite the above-mentioned symmetries, there is, to date, no known associated linear operator pair for 



the symmetric system. Without such a pair, the system (31a- 31b) cannot be solved by 1ST. The question 
remains, however, whether the symmetric system is integrable. In the absence of 1ST, the existence of 
multi-soliton solutions provides strong circumstantial evidence of integrability. 

In order to determine whether there are solitary waves which interact elastically- i.e. solitons, we first 



identify a solitary- wave solution of the symmetric system: a solution of the form (32) where the scalar q n is a 
one-soliton solution of IDNLS (as in eq. (^J)) is a solitary- wave solution of the symmetric system. Note that 
such a solution can have any polarization c subject only to the constraint that ||c|| 2 = 1. In order determine 
whether these vector solitary waves interact exactly, we use both direct and numerical methods. 

4.1 Soliton solutions by Hirota's Method 

Soliton solutions of the symmetric system can be found by Hirota's method. For concreteness, we consider 
the case where q„ has two components (N — 2). Under the independent variable transformation, 

gn 

q« = -j- 

In 

where g is a vector and /„ is a scalar, solutions of the bilinear equations 

ih 2 D t f n ■ g„ = /„_ig n+ i - 2/ n g„ + / n+ ig n _x (33a) 

fn+lfn-l- ti = h 2 \\g n \\ 2 (33b) 

are solutions of the symmetric system. 

The following solution of the bilinear equations yields a two-soliton solution: 

/ = 1 + e" 1 '"^.™ + e" 2 '"^." + |B 1 | 2 e " 1 .» + <» + " 2 -" + <™ (34a) 
g W = I (ghpi _ e -h P ^ e m, n (i + Bie V2, n +vl^ ( 34b ) 



1 (gh P 2 _ g-fcrf^ e r,2,n + 5 2e »l,n+<^ ( 34c ) 

where r)j, n and pj are as in eq. (|l^) and 
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Bi 



(e hpi - e hp2 ){e hpi + e hp z) 
( e h(pi+p 2 ) + i)( e /l (pi+P2) - 1) ' 



Bo = - 



(e 



hpx 



)(e hp i 



hp 2 



(e h( - p 



1+P2) 



In order to see that the above solution indeed gives a two-soliton solution, consider the long-time limits. 
We calculate the long-time limits of the Hirota form of the solution (34a -34c) by the same approa ch as wa s 
used for the asymmetric case. In the forward (+) and backward (— ) long-time limits, the solution a ( 34a -34c) 
asymptotically approaches the form 



as t — » ±00. 

The long-time limits for soliton 1 are 



<t>i \ sinh(aife) ipin 



J h 
+ _ (e*t\ wahjaih) „ 01 n 



where 







f>i = hbi, 



1 ' t '"sech(ai i „) 
e^sech^+loglSil) 



(35a) 
(35b) 



b\ = hb\ + arg B\ 



and ai,n, Pi,n are as m (17a- 17q ). In both limits, this is a solution of the form = c^qi i7l where q^ n is 
a one-soliton solution of IDNLS (as in (j|)). Note that the polarization, cf , is such that c^ 1 ^) 2 = 1 (and, 



therefore, 



,(2)±,2 



since lie 



±112 



1 



The long-time limits for soliton 2 are 



^t^e*."sech(a 2 ,„ + log| Bl |) 
h 



^ sinh(a 2 /i) j . 



e ip2 "sech(a 2 ,n) 



(36a) 
(36b) 



where 



hb 2 + arg B 2 



4 



hbo 



and a 2jn , $2,n are as in (17a 17b) with a 2 and & 2 replacing a\ and 61 respectively. In both the forward and 
backward long-time limit, this is a solution of the form c^q^ n where Ic^^l 2 = 1 (and, therefore, Ic^^l 2 = 



since c 



±112 



II 2 = !)• 



The solution (35a- 35b| . |36a -36b) is not the most general two-soliton interaction because 



0. 



(37) 



The two-soliton solution (35a -35b, 36a -36b) can be multiplied by a unitary matrix to obtain a two-soliton 
solution with any cj~,c^" such that |c^~ -c^"| =0 (under any such transformation, the condition |c+ • c£\ = 
will still hold). Thus, these two-soliton solutions are a constrained class of two-soliton solutions where ( |37| ) 
is the constraint. 
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More generally, Hirota's method can be used to derive solutions of the symmetric system with more 
than two solitons. These solutions can be represented as combinations of Pfaffians (these formulae are quite 
technical; they are presented in a separate publication [ fl4| ) . The M-soliton solutions derived in this manner 
are constrained such that either 



or 



c ■ • c fe | = 1 (38a) 



| C T. C fc| = (38b) 

(38c) 

for all j, k = 1, . . . , M. Moreover, the above satisfy 



(39) 



for all j = 1, . . . , M. These multi-soliton solutions are reminiscent of the special case of soliton interactions 
in the PDE discussed at the end of section 2.2: conditions ([38a -38b) and eq. ( |39| ) for the discrete symmetric 
system are the counterparts of (^a-9b) and eq. ( |l0| ) in the PDE. To date, there is no known analytic formula 



for a general multi-soliton solution of the symmetric system such that < c - ■ c k \ < 1 and |c ■ clj < 1. 



4.2 Numerical Simulation of Soliton Interactions 

The individual solitary waves of the symmetric system may have any polarization. In the PDE (^) such 
solitary waves interact as solitons. That is, there are vector multi-soliton solutions of VNLS with asymptotic 
polarizations such that < \cj ■ c^\ < 1 and \cj ■ c+| < 1. In the absence of analytical formulae in 
the discrete symmetric system for these general multi-soliton interactions, we investigated the collision of 
solitary waves by numerical simulation. 

In the simulations, the initial conditions were taken to be of the form 



q« =qi 



q 2 .« 



where q~ n = c J<jJ n an< ^ ijn 1S °f the f° rm ©• Then, the symmetric system was integrated in time by an 
adaptive Rungc-Kutta-Merson routine (from the NAG library) until the peaks were again well-separated (see 
Figure [I] for an example) . The separation of peaks in the initial and final conditions makes these conditions 
comparable to, respectively, the backwards and forwards long-time limits. The solitary-wave interactions 
simulated in this manner comprised initial conditions in which < |cj ■ c^7 1 < 1- Visually, in Figure |l| (and 
in other simulations) the solitary waves appear to interact without any radiation. We confirmed this finding 
by quantitative comparison of the solitary waves at the final time and exact solitary waves with the same 
height and speed parameters as the initial data. 

The error at the final time was defined separately for each solitary wave by 



Aj = —— max 



(40) 



where: 



q 



is the numerical data at the final time; 



q^ n is a solitary wave with the same amplitude and speed parameters- an and bj -as the initial data, q 



3,n> 



A 



sinh(ajh) 



is the amplitude of the exact solitary wave with amplitude parameter a,- ; 
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flj is the set of points containing the main contribution of the j-th solitary wave- i.e. flj = |n : | |q„ ^ 1 1 > e> 
for e small compared to max J= i,2 Aj. 

The resulting errors were small (see Table |l| for an example). Furthermore, when the user-supplied error 
bound in the adaptive integration scheme was decreased, the errors, Aj, decreased proportionally. Therefore, 
the differences between the final wave forms obtained by numerical simulation and exact solitary waves are 
accounted for by errors in the time integration. These results, strongly suggest that the solitary waves 
interact elastically- i. e. the solitary waves are solitons. 

In the numerical simulations described in Table [l], and in other experiments we considered initial data 
such that < |c^ • < 1, the general case of the soliton interaction for which there is no known analytical 
solution. In the PDE (||), such conditions result in the shift of the polarizations of the individual vector 
solitons- i.e. \cj ■ c+| < 1, j = 1, 2. The simulations described in Table [l[ consistent with other numerical 
experiments, show the same distinctive vector soliton behavior for the discrete symmetric system. 

The mechanism of the more general elastic soliton interactions observed in the symmetric system remains 
to be explained analytically. More generally, proof that the symmetric system is integrable remains as an 
important open problem. 
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Figure 1: Two-soliton interaction for the symmetric system. The filled boxes are \qn^\ and the open boxes 
are \qk |. Increasing time is read down column- wise. Soliton 1 is on the left and Soliton 2 is on the right at 
t = —8.46. They are reversed at t — 6.54. The soliton parameters are: a\ — 1, a 2 = 3, bi = .1, b 2 = — .1. 
The polarization vectors are: cj" = (1,0), = (.60, .80) , at t = -8.46, and c+ = (.33e i 237r , .95e l l27r ), 



c+ = (.84e- 



.54e" 



at t 



.204 < 1 and |c, • 



6.54. This is a typical two-soliton interaction where < |c 1 
= .993 < 1. 



.6 < 1, 
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a 2 = 1 


C 2 


Icr-cj-] 


l C 2 -CJI 


logio A i 


logio A 2 


(0.0,1.0) 


1.000 


1.000 


-7.57 


-7.57 


(0.2, .98) 


0.575 


0.575 


-7.54 


-7.54 


(0.4, .92) 


0.492 


0.492 


-7.52 


-7.52 


(0.6,0.8) 


0.624 


0.624 


-7.58 


-7.58 


(0.8,0.6) 


0.806 


0.806 


-7.56 


-7.56 


(1.0,0.0) 


1.000 


1.000 


-7.57 


-7.57 


0,2 = 3 


C 2 


|cr-cri 


l c 2 -cJI 


logio Ai 


logi A 2 


(0.0,1.0) 


1.000 


1.000 


-8.92 


-6.21 


(0.2, .98) 


0.900 


0.999 


-7.75 


-6.12 


(0.4, .92) 


0.617 


0.996 


-7.41 


-6.06 


(0.6,0.8) 


0.204 


0.993 


-7.27 


-6.07 


(0.8,0.6) 


0.396 


0.994 


-7.37 


-6.12 


(1.0,0.0) 


1.000 


1.000 


-7.61 


-6.25 



Table 1: Two-Soliton interaction for the symmetric system (Q). The soliton amplitude parameters are 
ai = 1 and a 2 — 1 or a 2 = 3 as it is given in each of the tables. For these values of a,j, the soliton width 
is comparable to the grid size and the solution is not "close" to the continuum limit. The soliton speed 
parameters, b\ = .1 and b 2 — — .1, are such that the solitons move slowly relative to one another thereby 
increasing the strength of the interaction. The polarizations before interaction are cj~ = (1,0) and c 2 as it 
is given in the tables. The polarizations after interaction are cf and cjj • Aj where j = 1,2, is the difference, 
as given by eq. (^), between the numerical solution at the final time and an exact solitary wave. The values 
\cj ■ Cj\ < 1 indicate that the polarization vectors are shifted by the soliton interaction. The small errors, 
Aj, show that the solitary waves interact nearly elastically with the measured deviation accounted for by 
error in the numerical time integration. 
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